#include "cxcore.h"
#include "cv.h"
#include "fftw3.h"

__declspec(dllexport) void RemovePixelBelowThreshold(IplImage* input, int threshold)
{
	int height, width, step, channels;
	uchar *src = 0;
	int i, j;

	// get the image data	
	height    = input->height;
	width     = input->width;
	step      = input->widthStep;
	channels  = input->nChannels;
	src      = (uchar*) input->imageData;

	for(i = 0; i < height; i++)
	{
		for(j = 0; j < width; j++)
		{
			if (src[i * step + j] <= threshold)
				src[i * step + j] = (uchar)0;
		}
	}
}

__declspec(dllexport) void Convert2Grayscale(IplImage* image, IplImage* output)
{

	//output = cvCreateImage(cvSize(image->width, image->height), 8, 1);
	cvCvtColor(image, output, CV_RGB2GRAY);
		
}

__declspec(dllexport) void InvertImage(IplImage* image, IplImage* img)
{
	//IplImage* img = cvCreateImage(cvSize(image->width, image->height), image->depth, image->nChannels);
	//cvCopyImage(image, img);

	int height, width, step, channels;
	uchar *data;
	int i, j, k;

	// get the image data	
	height    = img->height;
	width     = img->width;
	step      = img->widthStep;
	channels  = img->nChannels;
	data      = (uchar*) img->imageData;

	// invert the image
	for(i = 0; i < height; i++)
	{
		for(j = 0; j < width; j++)
		{
			for(k = 0; k < channels; k++)
			{
				data[i * step + j * channels + k] = 255 - data[i * step + j * channels + k];
			}
		}
	}

	//return img;
}


__declspec(dllexport) IplImage* BinarizeImage(IplImage* image, IplImage* )
{	
	if (image != NULL && image->nChannels == 1)
	{
		IplImage* resultImage = cvCreateImage(cvSize(image->width, image->height), 8, 1);
		cvThreshold(image, resultImage, 127, 255, cv::THRESH_OTSU);

		return resultImage;
	}
	else
	{
		return image;
	}
}

__declspec(dllexport) IplImage* ErodeImage(IplImage* image)
{
	IplImage* resultImage = cvCreateImage(cvSize(image->width, image->height), image->depth, image->nChannels);
	cvCopyImage(image,resultImage);

	cvErode(image, resultImage);

	return resultImage;
}

__declspec(dllexport) void SobelFilter(IplImage* image, IplImage* result) 
{ 
	/*IplImage* result = cvCreateImage(cvSize(image->width, image->height), image->depth, image->nChannels);
	cvCopyImage(image,result);*/

	if (image->nChannels == 1)
	{		
		cvFillImage(result, 0);

		int width = image->width;
		int height = image->height;

		int i, j;
		uchar* srcData = (uchar *)image->imageData;
		uchar* dstData = (uchar *)result->imageData;
		int step = image->widthStep/sizeof(uchar);

		int sumX, sumY, sum;

		for(j = 1; j < width - 2; j++) 
		{ 
			for(i = 1; i < height - 2; i++) 
			{ 
				sumX = -srcData[(i - 1) * step + j - 1] - 2 * srcData[(i - 1) * step + j] - srcData[(i - 1) * step + j + 1];
				sumX += srcData[(i + 1) * step + j - 1] + 2 * srcData[(i + 1) * step + j] + srcData[(i + 1) * step + j + 1];

				sumY = -srcData[(i - 1) * step + j - 1] - 2 * srcData[i * step + j - 1] - srcData[(i + 1) * step + j - 1];
				sumY += srcData[(i - 1) * step + j + 1] + 2 * srcData[i * step + j + 1] + srcData[(i + 1) * step + j + 1];

				sum = abs(sumX) + abs(sumY);
				if (sum > 255)
					sum = 255;
				else if (sum <= 30)
					sum = 0;
				dstData[i * step + j] = (uchar)sum;
			}	
		}
		//return result;
	}
	//return NULL;
}

__declspec(dllexport) IplImage* ResizeWithRate(IplImage* image, double rate)
{
	if (rate > 0.0)
	{
		int width = (int)(rate * image->width);
		int height = (int)(rate * image->height);
		IplImage *dest = cvCreateImage(cvSize(width, height), image->depth, image->nChannels);
		cvResize(image, dest);
		
		return dest;
	}

	return image;
}

__declspec(dllexport) IplImage* Resize(IplImage* image)
{
	double percent = 0.0;
	if (image->width > 1024)
	{
		percent = 1024.0 / (image->width * 1.0);
	}
	else if (image->height > 1024)
	{
		percent = 1024.0 / (image->height * 1.0);
	}

	if (percent != 0.0)
	{
		int width = (int)(percent * image->width);
		int height = (int)(percent * image->height);
		IplImage *dest = cvCreateImage(cvSize(width, height), image->depth, image->nChannels);
		cvResize(image, dest);
		
		return dest;
	}

	return image;
}

__declspec(dllexport) void EraseBackgroundWithThreshold1(IplImage* image, IplImage* result, int threshold)
{
	// Phong Le's code
	if (image->nChannels == 1)
	{
		cvCopyImage(image, result);

		// temp variables
		uchar* srcData = (uchar*)(image->imageData);
		uchar* desData = (uchar*)(result->imageData);

		int grayWidthStep = image->widthStep;

		// just keep pixel with value < threshold
		for (int y = 0; y < image->height; y++)
		{
			for (int x = 0; x < image->width; x++)
			{
				// erase pixel > threshold
				if (srcData[y * grayWidthStep + x] > threshold)
				{
					desData[y * grayWidthStep + x] = 0;
				}
			}
		}
	}
}

__declspec(dllexport) void EraseBackgroundWithThreshold2(IplImage* grayImage, IplImage* result, int backgray, int difThres)
{
	if (grayImage->nChannels != -1)
	{
		cvCopyImage(grayImage, result);

		// temp variables
		uchar* grayData = (uchar*)(result->imageData);

		int widthStep = result->widthStep;

		// just keep pixel with value < threshold
		for (int y = 0; y < grayImage->height; y++)
		{
			for (int x=0; x < grayImage->width; x++)
			{
				int i = y * widthStep + x;
				if (abs(grayData[i] - backgray) < difThres)
				{
					grayData[i] = 0;
				}
			}
		}
	}
}

__declspec(dllexport) void EraseBackground(IplImage* grayImage, IplImage* result)
{
	if (grayImage->nChannels != -1)
	{
		//IplImage* result  = cvCreateImage(cvSize(grayImage->width, grayImage->height), grayImage->depth, grayImage->nChannels);
		cvCopyImage(grayImage, result);

		// temp variables
		uchar* grayData = (uchar*)(result->imageData);

		int widthStep = result->widthStep;

		int grayscale [256];
		for (int i = 0; i < 256; i++)
		{
			grayscale[i] = 0;
		}
		//
		for (int y = 0; y < result->height; y++)
		{
			for (int x = 0; x < result->width; x++)
			{
				int i = y * widthStep + x;
				grayscale[grayData[i]]++;
			}
		}
		//
		int maxGray = 0;
		for (int i = 0; i < 256; i++)
		{
			if (grayscale[i] > grayscale[maxGray])
			{
				maxGray = i;
			}
		}
		// just keep pixel with value < threshold
		for (int y = 0; y < grayImage->height; y++)
		{
			for (int x=0; x < grayImage->width; x++)
			{
				int i = y * widthStep + x;
				if (abs(grayData[i] - maxGray) < 25)
				{
					grayData[i] = 0;
				}
			}
		}
	}
}

__declspec(dllexport) void RemoveGrid(IplImage* binaryImage, IplImage* resultImage, IplImage* gridImage)
{
	//IplImage* resultImage = cvCreateImage(cvSize(binaryImage->width, binaryImage->height), binaryImage->depth, binaryImage->nChannels);
	cvCopyImage(binaryImage, resultImage);

	cvSub(resultImage, gridImage, resultImage);
	cvSmooth(resultImage,resultImage, CV_MEDIAN,3);
	cvDilate(resultImage, resultImage, NULL, 1);

}

void removeFadeDetail(IplImage* input, int thresh)
{
	int height, width, step, channels;
	uchar *src = 0;
	int i, j;

	// get the image data	
	height    = input->height;
	width     = input->width;
	step      = input->widthStep;
	channels  = input->nChannels;
	src      = (uchar*) input->imageData;

	for(i = 0; i < height; i++)
	{
		for(j = 0; j < width; j++)
		{
			if (src[i * step + j] <= thresh)
				src[i * step + j] = (uchar)0;
		}
	}
}

__declspec(dllexport) void RemoveGrid_Fourier(IplImage* edgeImage, IplImage* resultImage, int filter_thickness)
{
	if (edgeImage->nChannels == 1)
	{
		uchar *edge_data;
		uchar *result_data;
		
		fftw_complex *data_in;
		fftw_complex *fft;        
		fftw_complex *ifft;    
		fftw_plan plan_f;
		fftw_plan plan_b;
	    
		int width, height, step;
		int i, j, k;

		/* get image properties */
		width = edgeImage->width;
		height = edgeImage->height;
		step = edgeImage->widthStep;
		edge_data = (uchar*)edgeImage->imageData;
		result_data = (uchar*)resultImage->imageData;
    
		/*initialize arrays for fftw operations */
		data_in = (fftw_complex*)fftw_malloc(sizeof(fftw_complex) * width * height );
		fft = (fftw_complex*)fftw_malloc(sizeof(fftw_complex) * width * height );
		ifft = (fftw_complex*)fftw_malloc(sizeof(fftw_complex) * width * height );
	
		/* create plans */
		plan_f = fftw_plan_dft_2d(height, width, data_in, fft, FFTW_FORWARD, FFTW_ESTIMATE);
		plan_b = fftw_plan_dft_2d(height, width, fft, ifft, FFTW_BACKWARD, FFTW_ESTIMATE);
   
		/* load edgeImage's data to fftw input */
		for(i = 0, k = 0 ; i < height ; i++)
		{
			for(j = 0 ; j < width ; j++)
			{
				data_in[k][0] = (double)edge_data[i * step + j];
				data_in[k][1] = 0.0;
				k++;
			}
		}
		/* perform FFT */
		fftw_execute(plan_f);
   		/* manipulate FTT */
		for(i = 0; i < filter_thickness; i++)
		{
			for(j = 0; j < width; j++)
			{
				fft[i * width + j][0] = (double)0;
				fft[i * width + j][1] = (double)0;

				fft[(height - 1 - i) * width + j][0] = (double)0;
				fft[(height - 1 - i) * width + j][1] = (double)0;
			}
		}
		/* perform IFFT */
		fftw_execute(plan_b);
	    
		/* normalize IFFT result */
		for(i = 0 ; i < (width * height) ; i++)
		{
			ifft[i][0] /= (double)(width * height);
		}
	    
		/* copy IFFT result to result data */
		for( i = 0, k = 0 ; i < height ; i++ ) 
		{
			for( j = 0 ; j < width ; j++ ) 
			{
				double t = ifft[k++][0];
				if (t < 0.0)
					t = 0.0;
				else if (t > 255)
					t = 255;
				result_data[i * step + j] = (uchar)t;
			}
		}
    
		/* free memory */
		fftw_destroy_plan(plan_f);
		fftw_destroy_plan(plan_b);
		fftw_free(data_in);
		fftw_free(fft);
		fftw_free(ifft);
		removeFadeDetail(resultImage, 30);
	}
}

__declspec(dllexport) void GetGrid(IplImage* binaryImage, IplImage* resultImage, double angle)
{
	if (binaryImage->nChannels == 1)
	{
		double cos_thres = cos(angle * CV_PI / 180);
		double module = 0.0;
		double horizontalCos = 0.0;
		double verticalCos = 0.0;
		
		// detect line using HoughLine transform
		CvMemStorage* storage2 = cvCreateMemStorage(0);
		CvSeq* lines = cvHoughLines2(binaryImage, storage2,
			CV_HOUGH_PROBABILISTIC, 1, CV_PI/180, 50, 50, 25 );

		//IplImage* resultImage = cvCreateImage(cvSize(binaryImage->width, binaryImage->height), binaryImage->depth, binaryImage->nChannels);
		cvFillImage(resultImage, 0);

		for(int i = 0; i < lines->total; i++ )
		{
			CvPoint* line = (CvPoint*)cvGetSeqElem(lines, i);
			module = sqrt(pow(1.0 * line[1].x - line[0].x, 2) + pow(1.0 * line[1].y - line[0].y, 2));
			horizontalCos = 1.0 * abs(line[1].x - line[0].x) / module;
			verticalCos = 1.0 * abs(line[1].y - line[0].y) / module;

			if (horizontalCos > cos_thres || verticalCos > cos_thres)
			{
				if (module /binaryImage->width >= 0.15 || module /binaryImage->height >= 0.15)
				{
					cvLine(resultImage, line[0], line[1], CV_RGB(255,255,255), 1, 8 );
				}
			}
		}
	}
}

void increaseBrightness(IplImage* input, int amount)
{
	int height, width, step, channels;
	uchar *src = 0;
	int i, j;

	// get the image data	
	height    = input->height;
	width     = input->width;
	step      = input->widthStep;
	channels  = input->nChannels;
	src      = (uchar*) input->imageData;

	for(i = 0; i < height; i++)
	{
		for(j = 0; j < width; j++)
		{
			src[i * step + j] += (uchar)amount;
		}
	}
}

__declspec(dllexport) void GetGrid_Fourier(IplImage* edgeImage, IplImage* resultImage, int filter_thickness)
{
	if (edgeImage->nChannels == 1)
	{
		uchar *edge_data;
		uchar *result_data;
		
		fftw_complex *data_in;
		fftw_complex *fft;        
		fftw_complex *ifft;    
		fftw_plan plan_f;
		fftw_plan plan_b;
	    
		int width, height, step;
		int i, j, k;

		/* get image properties */
		width = edgeImage->width;
		height = edgeImage->height;
		step = edgeImage->widthStep;
		edge_data = (uchar*)edgeImage->imageData;
		result_data = (uchar*)resultImage->imageData;
    
		/*initialize arrays for fftw operations */
		data_in = (fftw_complex*)fftw_malloc(sizeof(fftw_complex) * width * height );
		fft = (fftw_complex*)fftw_malloc(sizeof(fftw_complex) * width * height );
		ifft = (fftw_complex*)fftw_malloc(sizeof(fftw_complex) * width * height );
	
		/* create plans */
		plan_f = fftw_plan_dft_2d(height, width, data_in, fft, FFTW_FORWARD, FFTW_ESTIMATE);
		plan_b = fftw_plan_dft_2d(height, width, fft, ifft, FFTW_BACKWARD, FFTW_ESTIMATE);
   
		/* load edgeImage's data to fftw input */
		for(i = 0, k = 0 ; i < height ; i++)
		{
			for(j = 0 ; j < width ; j++)
			{
				data_in[k][0] = (double)edge_data[i * step + j];
				data_in[k][1] = 0.0;
				k++;
			}
		}
		/* perform FFT */
		fftw_execute(plan_f);
   		/* manipulate FTT */
		for(i = 0; i < filter_thickness; i++)
		{
			for(j = 1; j < width - 1; j++)
			{
				fft[i * width + j][0] = (double)0;
				fft[i * width + j][1] = (double)0;

				fft[(height - 1 - i) * width + j][0] = 0.0;
				fft[(height - 1 - i) * width + j][1] = 0.0;
			}
		}
		/* perform IFFT */
		fftw_execute(plan_b);
	    
		/* normalize IFFT result */
		for(i = 0 ; i < (width * height) ; i++)
		{
			ifft[i][0] /= (double)(width * height);
		}
	    
		/* copy IFFT result to result data */
		for(i = 0, k = 0 ; i < height ; i++) 
		{
			for(j = 0 ; j < width ; j++) 
			{
				result_data[i * step + j] = (uchar)ifft[k++][0];
			}
		}
    
		/* Get grid */	
		cvSub(edgeImage, resultImage, resultImage);
		increaseBrightness(resultImage, 15);	
		/* free memory */
		fftw_destroy_plan(plan_f);
		fftw_destroy_plan(plan_b);
		fftw_free(data_in);
		fftw_free(fft);
		fftw_free(ifft);
	}
}

__declspec(dllexport) void ContrastStretch(IplImage* grayImage, IplImage* resultImage)
{
	//IplImage* resultImage = cvCreateImage(cvSize(grayImage->width, grayImage->height), grayImage->depth, grayImage->nChannels);
	//cvCopyImage(grayImage, resultImage);

	cvNormalize(grayImage, resultImage, 0, 255, CV_MINMAX); 

}

__declspec(dllexport) void BlurGrid(IplImage* grayImage, IplImage* resultImage, IplImage* binaryGridImage)
{
	if (grayImage->nChannels == 1 && binaryGridImage->nChannels == 1)
	{
		//IplImage* resultImage = cvCreateImage(cvSize(grayImage->width, grayImage->height), grayImage->depth, grayImage->nChannels);
		cvCopyImage(grayImage, resultImage);

		// temp variables
		uchar* grayData = (uchar*)(resultImage->imageData);
		uchar* binaryData = (uchar*)(binaryGridImage->imageData);

		int widthStep = resultImage->widthStep;
		// 
		int grayscale [256];
		for (int i = 0; i < 256; i++)
		{
			grayscale[i] = 0;
		}
		//
		for (int y = 0; y < resultImage->height; y++)
		{
			for (int x = 0; x < resultImage->width; x++)
			{
				int i = y * widthStep + x;
				grayscale[grayData[i]]++;
			}
		}
		// 
		int maxGray = 0;
		for (int i = 0; i < 256; i++)
		{
			if (grayscale[i] > grayscale[maxGray])
			{
				maxGray = i;
			}
		}		
		//
		for (int y = 0; y < resultImage->height; y++)
		{
			for (int x = 0; x < resultImage->width; x++)
			{
				int i = y * widthStep + x;
				if (binaryData[i] == 255)
				{
					if (abs(grayData[i] - maxGray) < 20)
					{
						grayData[i] = maxGray;
					}
				}
			}
		}
	}
}